#=
测试高斯积分
=#


include("../src/lgvmary.jl")
using ..lgvmary

using Symbolics


#=
∫ e^{-x^2} dx = √π
∫ x^{2n} e^{-0.5ax^2} dx = (2n-1)!! sqrt(2π/a) / a^n 
n=1 a=2: √π/2
∫ x^{2} e^{-x^2} dx / ∫ e^{-x^2} dx = 1/2
#
integral from -inf to inf x^2 e^{-x^2 + i*x} dx = √π/(4*e^{1/4}) ≈ 0.345097
int e^{-x^2 + i*x}  from -inf to inf = √π/e^{1/4}
∫ x^{2} e^{-x^2+x*im} dx / ∫ e^{-x^2+x*im} dx = 1/4
=#


#@variables x[1:3]

#println(x, x[1], x[2])

function run()
    x, imagv, S = dissipative_action_Gaussian()
    O = x*exp(S)
    ∂S_∂x = Symbolics.derivative(S, x)
    #expand_derivatives(Differential(x)(S))

    println("$x $S $∂S_∂x $O")

    Sval = Symbolics.value(substitute(S, Dict(x => 5, imagv => 0.0+1.0im)))
    #Symbolics.fixpoint_sub(S, Dict(x => 5))
    Spval = Symbolics.value(substitute(∂S_∂x, Dict(x => 5, imagv => 0.0+1.0im)))
    #Spval = Symbolics.fixpoint_sub(∂S_∂x, Dict(x => 5, imagv => 0.0+1.0im))

    println("$Sval $Spval")


    xmin = -10
    xmax = 10
    sct = 100000
    sps = (xmax-xmin)*rand(sct) .+ xmin
    num = 0.0
    den = 0.0
    eS = exp(-S)
    eSO = eS*x^2
    println("$eS $eSO")
    for v in sps
        num += Symbolics.value(Symbolics.substitute(eSO, Dict(x=>v, imagv => 0.0+1.0im)))
        den += Symbolics.value(Symbolics.substitute(eS, Dict(x=>v, imagv => 0.0+1.0im)))
    end
    println(num/den)
end


run()
